Added method for gradientSymTensor by qiUip · Pull Request #127 · ucl-bug/k-wave-ii

qiUip · 2024-06-17T14:35:52Z

See issue #124 for description and #120 for the solver that requires this method.

The PR is made with awareness that the API might have to change. It is also worth exploring the further use of the gradient of a symmetric tensor operator as I'm not sure if it will be required for any other solvers.

github-actions · 2024-06-17T14:38:29Z

Cyclomatic Complexity Check

L 97 (C 22-41): The McCabe cyclomatic complexity of 'initialiseSensorData' is 19.
L 153 (C 24-39): The McCabe cyclomatic complexity of 'recordSensorData' is 6.
L 159 (C 24-37): The McCabe cyclomatic complexity of 'AcousticSolver' is 5.
L 189 (C 29-40): The McCabe cyclomatic complexity of 'get.pressure' is 1.
L 192 (C 33-48): The McCabe cyclomatic complexity of 'get.densitySplit' is 1.
L 195 (C 29-40): The McCabe cyclomatic complexity of 'get.velocity' is 1.
L 38 (C 21-39): The McCabe cyclomatic complexity of 'autoComputeTimeStep' is 3.
L 34 (C 10-24): The McCabe cyclomatic complexity of 'executeTimeStep' is 25.
L 60 (C 7): The McCabe cyclomatic complexity of the anonymous function on line 60 is 1.
L 61 (C 9): The McCabe cyclomatic complexity of the anonymous function on line 61 is 1.
L 62 (C 12): The McCabe cyclomatic complexity of the anonymous function on line 62 is 1.
L 63 (C 14): The McCabe cyclomatic complexity of the anonymous function on line 63 is 1.
L 67 (C 54-77): The McCabe cyclomatic complexity of 'buildAbsorptionOperators' is 2.
L 33 (C 23-42): The McCabe cyclomatic complexity of 'kappaSplitCorrection' is 2.
L 46 (C 10-34): The McCabe cyclomatic complexity of 'setAbsorptionCoefficients' is 5.

L 59 (C 14-17): The McCabe cyclomatic complexity of 'curl' is 2.
L 64 (C 14-23): The McCabe cyclomatic complexity of 'divergence' is 1.
L 66 (C 14-28): The McCabe cyclomatic complexity of 'divergenceSplit' is 16.
L 72 (C 15-35): The McCabe cyclomatic complexity of 'divergenceTensorSplit' is 20.
L 58 (C 16-28): The McCabe cyclomatic complexity of 'fracLaplacian' is 4.
L 64 (C 15-22): The McCabe cyclomatic complexity of 'gradient' is 15.
L 60 (C 17-23): The McCabe cyclomatic complexity of 'stagger' is 17.
L 52 (C 10-21): The McCabe cyclomatic complexity of 'validateSize' is 7.
L 68 (C 24-43): The McCabe cyclomatic complexity of 'validateNonNegScalar' is 2.
L 84 (C 24-41): The McCabe cyclomatic complexity of 'validateIndexField' is 6.
L 99 (C 24-32): The McCabe cyclomatic complexity of 'Materials' is 1.
L 124 (C 32-42): The McCabe cyclomatic complexity of 'addMaterial' is 20.
L 222 (C 22-34): The McCabe cyclomatic complexity of 'listMaterials' is 4.
L 254 (C 22-40): The McCabe cyclomatic complexity of 'listMaterialIndices' is 1.
L 263 (C 26-41): The McCabe cyclomatic complexity of 'getMaterialNames' is 1.
L 270 (C 45-76): The McCabe cyclomatic complexity of 'getMaterialNamesAndIndicesSorted' is 2.
L 128 (C 24-29): The McCabe cyclomatic complexity of 'Medium' is 1.
L 141 (C 24-31): The McCabe cyclomatic complexity of 'subsasgn' is 5.
L 171 (C 24-30): The McCabe cyclomatic complexity of 'refresh' is 1.
L 179 (C 18-40): The McCabe cyclomatic complexity of 'set.soundSpeedReference' is 2.
L 185 (C 18-39): The McCabe cyclomatic complexity of 'set.diffusionReference' is 2.
L 196 (C 24-46): The McCabe cyclomatic complexity of 'updateDerivedPaddedMaps' is 10.
L 246 (C 25-35): The McCabe cyclomatic complexity of 'mapProperty' is 6.
L 284 (C 22-40): The McCabe cyclomatic complexity of 'get.absorptionPower' is 4.
L 57 (C 36-56): The McCabe cyclomatic complexity of 'get.totalSensorPoints' is 1.
L 67 (C 33-49): The McCabe cyclomatic complexity of 'ProcessSensorData' is 4.
L 87 (C 24-29): The McCabe cyclomatic complexity of 'Solver' is 5.
L 81 (C 24-39): The McCabe cyclomatic complexity of 'TimeDomainSolver' is 1.
L 110 (C 18-37): The McCabe cyclomatic complexity of 'setInitialConditions' is 1.
L 57 (C 10-12): The McCabe cyclomatic complexity of 'run' is 5.
L 1 (C 1): The McCabe cyclomatic complexity is 0.

codecov · 2024-06-17T14:39:20Z

Codecov Report

❌ Patch coverage is 76.00000% with 18 lines in your changes missing coverage. Please review.
✅ Project coverage is 98.36%. Comparing base (c474ac2) to head (491ebcd).
⚠️ Report is 28 commits behind head on develop.

Files with missing lines	Patch %	Lines
...oolbox/@FourierCollocation/divergenceTensorSplit.m	76.00%	18 Missing ⚠️

Additional details and impacted files

@@             Coverage Diff             @@
##           develop     #127      +/-   ##
===========================================
- Coverage   100.00%   98.36%   -1.64%     
===========================================
  Files           59       48      -11     
  Lines          911     1102     +191     
===========================================
+ Hits           911     1084     +173     
- Misses           0       18      +18

Flag	Coverage Δ
unit	`98.36% <76.00%> (-1.64%)`	⬇️

Flags with carried forward coverage won't be shown. Click here to find out more.

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:

❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

btreeby

@qiUip Ben and I spent some time working through the code and assumptions. Broadly looks good! See comments about name and input ordering, etc.

btreeby · 2024-06-19T08:59:21Z

+% components are [sxx] in 1D, [sxx, syy, sxy]' in 2D and 
+% [sxx, syy, szz, sxy, sxz, syz]' in 3D.


We should use the standard Voigt notation order here unless there is a good reason not to, which is [s11, s22, s33, s23, s13, s12] or [sxx, syy, szz, syz, sxz, sxy]. (This will change the indexing used later on.)

Thanks for the suggestion! The reason I did not adopt the Voigt notation was that it made it less intuitive (at least as an initial thought) to apply the Einstein notation in a loop in the solver (see this line).

In the loop, using the indices for i = 1 $(x)$ and j = 2 $(y)$ of the velocity gradient tensor, adding 1 makes it index 4 easily for the $\sigma_{xy}$ component. I'm sure there is a way to adapt the indexing in the constitutive equation summation loops to use the Voigt notation, but I will have to have a bit of a think how to best go about it.

btreeby · 2024-06-19T09:10:25Z

@@ -0,0 +1,176 @@
+%% gradientSymTensor


The operation we're applying here is actually the divergence of a tensor, but in this case, we're keeping the split components. So I'd call this divergenceSplitTensor (similar in spirit to divergenceSplit).

For completeness, we can also write a divergenceTensor which just calls divergenceSplitTensor and then sums the components over the 5th dimension, leaving a vector field (see also divergence).

I would leave out the sym part in the function name. Later we could add support for non-symmetric tensors, where the input is given as a Tensor (rather than a vector with the components in Voigt notation). The output will be the same size.

The docs will need updating (where you write gradient).

Sounds good! Thanks for that. I will address the comments and the name change at the same time next chance I get.

I think what threw me off here is equation (7a) in the IEEE paper and I just didn't pay enough attention to the maths.

(1) Calculate the spatial gradients of the stress field using the Fourier collocation spectral method.

I have applied this change. Please carefully read through the description as I may have gotten the terminology wrong.

Also note that I didn't (yet) add the methods that will sum-up the components to a vector. I will get around to that in the near future, but first I wanted to make sure what's here is correct. I'm not sure if we should use the split directly in the solver or use the summed version due to the multi-axisal pml - I will look into that as soon as I can.

btreeby · 2024-06-19T12:05:59Z

+            case 3
+                % Indecies for the symmetric (squashed) stress tensor are
+                % different in 2D and 3D so the spatial derivatives for the
+                % dxy components need overwritten with the correct


I would handle the case statement so that you don't have to compute these twice (the FFTs are expensive and this operation is done many times). Same for below for case with ND k-space correction.

Very good point! I'll address this by the weekend as soon as I find some time.

I used an if statement inside the case 2 statement. This is a simple solution but I think it would work as it is outside the internal loops over the grid points.

btreeby · 2024-06-19T12:53:25Z

+        dsxxdx = obj.ddxShiftPos;
+        dsxydy = obj.ddyShiftNeg;
+        dsxzdz = obj.ddzShiftNeg;
+        dsxydx = obj.ddxShiftNeg;
+        dsyydy = obj.ddyShiftPos;
+        dsyzdz = obj.ddzShiftNeg;
+        dsxzdx = obj.ddxShiftNeg;
+        dsyzdy = obj.ddyShiftNeg;
+        dszzdz = obj.ddzShiftPos;


Noting that the off-diagonal components are duplicated, but I can't think of a nice notation to combine them, and it shouldn't have any performance or memory consequences anyway.

Do you mean that obj.ddxShiftNeg is called twice, once for dsxydx and once for dsxzdx (and the same for the y and z components)? The other option is to set dsxzdx = dsxydx then, but I'm not sure I like that either.

…ed fft calculations in 3D

#124: Added method for gradientSymTensor with unit tests

27c9d11

btreeby requested changes Jun 19, 2024

View reviewed changes

Mashy Green and others added 4 commits June 19, 2024 16:08

#124: Correcting some typos in the comments

5a5fdd0

#124: Changed the name to diverganceTensorSplit and removed duplicatt…

ea52253

…ed fft calculations in 3D

#124: Updated unit tests to account for name change

cd8a12f

Merge remote-tracking branch 'origin/develop' into 124-gradientSymTensor

491ebcd

		% components are [sxx] in 1D, [sxx, syy, sxy]' in 2D and
		% [sxx, syy, szz, sxy, sxz, syz]' in 3D.

Uh oh!

Conversation

qiUip commented Jun 17, 2024

Uh oh!

github-actions Bot commented Jun 17, 2024 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Cyclomatic Complexity Check

Uh oh!

codecov Bot commented Jun 17, 2024 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Codecov Report

Uh oh!

btreeby left a comment

Choose a reason for hiding this comment

Uh oh!

Uh oh!

btreeby Jun 19, 2024

Choose a reason for hiding this comment

Uh oh!

qiUip Jun 19, 2024 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

btreeby Jun 19, 2024

Choose a reason for hiding this comment

Uh oh!

qiUip Jun 19, 2024 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

qiUip Jun 19, 2024

Choose a reason for hiding this comment

Uh oh!

Uh oh!

Uh oh!

btreeby Jun 19, 2024

Choose a reason for hiding this comment

Uh oh!

qiUip Jun 19, 2024

Choose a reason for hiding this comment

Uh oh!

qiUip Jun 19, 2024

Choose a reason for hiding this comment

Uh oh!

btreeby Jun 19, 2024

Choose a reason for hiding this comment

Uh oh!

qiUip Jun 19, 2024 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

3 participants

github-actions Bot commented Jun 17, 2024 •

edited

Loading

codecov Bot commented Jun 17, 2024 •

edited

Loading

qiUip Jun 19, 2024 •

edited

Loading

qiUip Jun 19, 2024 •

edited

Loading

qiUip Jun 19, 2024 •

edited

Loading